#include <iostream>
#include "LHEF.h"
#include <iomanip>
#include <fstream>
#include <iterator>
#include <vector>
#include <sstream>
#include <string>
#include <string.h>
#include <cstdlib>
#include <stdio.h>
#include <stdio.h>
#include <stdlib.h>
#include <math.h>

using namespace std;

//*************************************************
vector<string> getCSstrings(string inputstring) {
	//*************************************************
	vector<string> result;
	istringstream instream;
	instream.str(inputstring);
	string tempstring;

	while (getline(instream, tempstring, ',')) {
		//cout << tempstring << endl;
		result.push_back(tempstring);
	}
	return result;
}

//*********************************************
void writeLHEeventline(ofstream & homemade_lhe, LHEF::HEPEUP hepeup, int npart) {
	//*********************************************

	homemade_lhe << npart << " " << hepeup.IDPRUP << " " << hepeup.XWGTUP << " " << hepeup.SCALUP << " " << hepeup.AQEDUP << " " << hepeup.AQCDUP << endl;
}

//*********************************************
void writeLHEeventlineRescaled(ofstream & homemade_lhe, LHEF::HEPEUP hepeup, int npart, float rescaling) {
	//*********************************************

	homemade_lhe << npart << " " << hepeup.IDPRUP << " " << hepeup.XWGTUP << " " << 1.0 * rescaling * hepeup.SCALUP << " " << hepeup.AQEDUP << " " << hepeup.AQCDUP << endl;
}

//*********************************************
void writeLHEeventlineReweight(ofstream & homemade_lhe, LHEF::HEPEUP hepeup, int npart, float reweighting) {
	//*********************************************
    
	homemade_lhe << npart << " " << hepeup.IDPRUP << " " << 1.0 * reweighting * hepeup.XWGTUP << " " <<  hepeup.SCALUP << " " << hepeup.AQEDUP << " " << hepeup.AQCDUP << endl;
}

//*********************************************
void writeLHEeventlineNewWeight(ofstream & homemade_lhe, LHEF::HEPEUP hepeup, int npart, float newweight) {
	//*********************************************
    
	homemade_lhe << npart << " " << hepeup.IDPRUP << " " << newweight << " " <<  hepeup.SCALUP << " " << hepeup.AQEDUP << " " << hepeup.AQCDUP << endl;
}

//*********************************************
void writeLHEeventlineNewScale(ofstream & homemade_lhe, LHEF::HEPEUP hepeup, int npart, float newscale) {
	//*********************************************

	homemade_lhe << npart << " " << hepeup.IDPRUP << " " << hepeup.XWGTUP << " " <<  newscale << " " << hepeup.AQEDUP << " " << hepeup.AQCDUP << endl;
}


//*********************************************
void writeLHEparticletoLHE_OriginalStatus(ofstream & homemade_lhe, LHEF::HEPEUP hepeup, int i) {
	//*********************************************
	homemade_lhe << hepeup.IDUP[i] << " " << hepeup.ISTUP[i] << " " << hepeup.MOTHUP[i].first << " " << hepeup.MOTHUP[i].second << " " << hepeup.ICOLUP[i].first << " " << hepeup.ICOLUP[i].second << " " << hepeup.PUP[i][0] << " " << hepeup.PUP[i][1] << " " << hepeup.PUP[i][2] << " " << hepeup.PUP[i][3] << " " << hepeup.PUP[i][4] << " " << hepeup.VTIMUP[i] << " 9" << endl;
}

//*********************************************
void writeLHEparticletoLHEChangePID(ofstream & homemade_lhe, LHEF::HEPEUP hepeup, int i, int pidin, int pidout) {
	//*********************************************
	int newpid;
	if (hepeup.IDUP[i] == pidin) {
		//hepeup.IDUP[i] = pidout;

		newpid = pidout;
//		cout<< " found a " << hepeup.IDUP[i] << "=="<<  pidin << "---> Changed to:"<< newpid <<endl;
	} else {
		newpid = hepeup.IDUP[i];
	}
	homemade_lhe << newpid << " " << hepeup.ISTUP[i] << " " << "0" << " " << "0" << " " << hepeup.ICOLUP[i].first << " " << hepeup.ICOLUP[i].second << " " << hepeup.PUP[i][0] << " " << hepeup.PUP[i][1] << " " << hepeup.PUP[i][2] << " " << hepeup.PUP[i][3] << " " << hepeup.PUP[i][4] << " " << hepeup.VTIMUP[i] << " 9" << endl;
}

//*********************************************
void writeLHEparticletoLHE(ofstream & homemade_lhe, LHEF::HEPEUP hepeup, int i) {
	//*********************************************
	homemade_lhe << hepeup.IDUP[i] << " " << "1" << " " << "0" << " " << "0" << " " << hepeup.ICOLUP[i].first << " " << hepeup.ICOLUP[i].second << " " << hepeup.PUP[i][0] << " " << hepeup.PUP[i][1] << " " << hepeup.PUP[i][2] << " " << hepeup.PUP[i][3] << " " << hepeup.PUP[i][4] << " " << hepeup.VTIMUP[i] << " 9" << endl;
}

int main(int argc, char ** argv) {

	char blue[] = { 0x1b, '[', '1', ';', '3', '4', 'm', 0 };
	char ublue[] = { 0x1b, '[', '4', ';', '3', '4', 'm', 0 };
	char rblue[] = { 0x1b, '[', '7', ';', '3', '4', 'm', 0 };
	char magenta[] = { 0x1b, '[', '1', ';', '3', '5', 'm', 0 };
	char red[] = { 0x1b, '[', '1', ';', '3', '1', 'm', 0 };
	char normal[] = { 0x1b, '[', '0', ';', '3', '9', 'm', 0 };

	//	cout<<argc;
	if (argc < 3) {
		puts("usage: lhetools merge/rescale/reweight/newscale/newweight/changepid/targetcrosssection inputfile [outputfile] ");
		puts("usage: lhetools merge inputfile1 inputfile2 ... inputfileN outputfile ");
		puts("usage: lhetools rescale,2 inputfile [outputfile] ");
        puts("usage: lhetools reweight,1.2 inputfile [outputfile] ");
		puts("usage: lhetools newscale,91.18 inputfile [outputfile] ");
        puts("usage: lhetools newweight,1.0 inputfile [outputfile] ");
		puts("usage: lhetools cut,1000 inputfile [outputfile] ");
		puts("usage: lhetools changepid,11,13 inputfile [outputfile] ");
        puts("usage: lhetools changepid,11,13 inputfile [outputfile] ");
		exit(1);
	}


    
	//***************************************************
	/*
	 double pTmin;
	 double dcutmin;
	 pTmin=20;
	 dcutmin=400;
	 */
	//***************************************************


	char* inputfile;
	char* outfile;
	char* inclexcl;
    

	string input_action;
	input_action = argv[1];
	vector<string> action = getCSstrings(input_action);
	string operation = action[0];
	float rescaling,newscale,reweighting,newweight,targetcrosssection,targetweigth;
	int pidin, pidout;
	int endpoint = 0;

	if (operation == "merge") {
		// do the merging
	} else if (operation == "rescale") {
		rescaling = atof(action[1].c_str());
	} else if (operation == "reweight") {
		reweighting = atof(action[1].c_str());
	} else if (operation == "newscale") {
		newscale = atof(action[1].c_str());
	} else if (operation == "newweight") {
		newweight = atof(action[1].c_str());
	} else if (operation == "targetcrosssection") {
		targetcrosssection = atof(action[1].c_str());
	}else if (operation == "changepid") {
		pidin = atoi(action[1].c_str());
		pidout = atoi(action[2].c_str());
	} else if (operation == "cut") {
		endpoint = atoi(action[1].c_str());
	} else {
		puts("no valid operation has been asked");
		exit(1);
	}

    if (argc == 3) {
        cout<<  "output file will be assigned by default" <<  endl;
        
        std::string str = argv[2] + input_action + ".lhe";
        
        char *cstr = new char[str.length() + 1];
        strcpy(cstr, str.c_str());
        
        outfile = cstr;
 
    } else
    {
        outfile = argv[argc - 1];
    }
    
    cout <<  "will output to "<< red << outfile << normal << endl;

    
    
	ofstream homemade_lhe;
	homemade_lhe.open(outfile);
   
	homemade_lhe << endl;
	homemade_lhe.close();
	LHEF::Writer writer(outfile);
    
	/*begin operation on the inputfiles */
	for (int file = 2; (file < argc - 1 && argc > 3 ) || (file < argc && argc == 3) ; file++) {
        puts("back here");
		inputfile = argv[file];
		LHEF::Reader countreader(inputfile);
        LHEF::Reader reader(inputfile);
        
        

		if (file == 2) {
			// Copy header and init blocks and write them out.
			//  if ( reader.outsideBlock.length() ) std::cerr << reader.outsideBlock;
			writer.headerBlock() << reader.headerBlock;
			writer.initComments() << reader.initComments;
			writer.heprup = reader.heprup;
			writer.init();

			homemade_lhe.open(outfile, ios::app);
			homemade_lhe << "                                             " << endl;
			homemade_lhe.close();
		}
        
		long numberofevents = 0;

        while (countreader.readEvent() == true && (numberofevents < endpoint || endpoint == 0)) {
            
			numberofevents++;
		}
        cout << "the sample contains " << red << numberofevents << normal << " events"<< endl;
        
        if (operation == "targetcrosssection") {
            targetweigth = targetcrosssection/(1.0 * numberofevents );
            cout << "events will be reweighted to a constant weight w="<< targetweigth <<" so that the entire sample will have total cross-section "<< targetcrosssection <<" pb"<< endl;
        }

        
        long neve = 0;
        
        
		while (reader.readEvent() == true && (neve < endpoint || endpoint == 0)) {

			if (reader.outsideBlock.length()) {
				//			std::cout << reader.outsideBlock;
			}
			//cout << neve << endl;

			homemade_lhe.open(outfile, ios::app);
			homemade_lhe << "<event>" << endl;
			//cout << "writing event description of event # " << neve << endl;

			if (operation == "rescale") {
				writeLHEeventlineRescaled(homemade_lhe, reader.hepeup, reader.hepeup.NUP, rescaling);
			} else if (operation == "newscale") {
				 writeLHEeventlineNewScale(homemade_lhe, reader.hepeup, reader.hepeup.NUP, newscale);
			} else if (operation == "newweight") {
                writeLHEeventlineNewWeight(homemade_lhe, reader.hepeup, reader.hepeup.NUP, newweight);
			} else if (operation == "targetcrosssection") {
                writeLHEeventlineNewWeight(homemade_lhe, reader.hepeup, reader.hepeup.NUP, targetweigth);
			} else if (operation == "reweight") {
                writeLHEeventlineReweight(homemade_lhe, reader.hepeup, reader.hepeup.NUP, reweighting);
			}
            
				else {
				writeLHEeventline(homemade_lhe, reader.hepeup, reader.hepeup.NUP);
			}
			for (int particle = 0; particle < reader.hepeup.NUP; particle++) {
                
                    // cout << particle << endl;

				// deal with paricles copy and call the approriate writer according to the called operation

				if (operation == "changepid") {
					writeLHEparticletoLHEChangePID(homemade_lhe, reader.hepeup, particle, pidin, pidout);
				} else if (operation=="cut" || operation == "merge" || operation == "rescale" || operation=="newscale" || operation == "reweight" || operation == "newweight" || operation == "targetcrosssection") {
					writeLHEparticletoLHE_OriginalStatus(homemade_lhe, reader.hepeup, particle);
				}
			}

			homemade_lhe << "</event>\n" << endl;
			homemade_lhe.close();

			neve++;
		}

		cout << inputfile << " has been processed" << endl;
	}
    
	homemade_lhe.open(outfile, ios::app);
    
	homemade_lhe << "</LesHouchesEvents>" << endl;
	homemade_lhe.close();
    cout << outfile << " has been closed" << endl;

	//	writer.LHEF::Writer::~Writer();
	//cout << extra_lepton << " events have at least an extra lepton due to showering" << endl;

}
